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Abstract In some studies requiring predictive and CPU- time consuming numerical 
models, the sampling design of the model input variables has to be chosen with cau- 
tion. For this purpose, Latin hypercube sampling has a long history and has shown its 
robustness capabilities. In this paper we propose and discuss a new algorithm to build 
a Latin hypercube sample (LHS) taking into account inequality constraints between 
the sampled variables. This technique, called constrained Latin hypercube sampling 
(cLHS), consists in doing permutations on an initial LHS to honor the desired mono- 
tonic constraints. The relevance of this approach is shown on a real example concern- 
ing the numerical welding simulation, where the inequality constraints are caused by 
the physical decreasing of some material properties in function of the temperature. 

Keywords Computer experiment • Latin hypercube sampling • Design of Experi- 
ments • Uncertainty analysis • Dependence 



1 Introduction 

With the advent of computing technology and numerical methods, investigation of 
computer code experiments remains an important challenge. Complex computer mod- 
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els calculate several output values which can depend on a large number of input 
parameters and physical variables. These computer models are used to make sim- 
ulations as well as p redictions, uncertainty and sensitivity analyses or to solve op 



timization problems QFang 
Lew and Steinber j, 2010l) . 



timization problems dFang et al.L l2006t iKleiinenl l2008t iDe Rocquignv et all 12008 



However, complex computer codes are often too time expensive to be directly 
used to perform such studies. For uncertainty propagation and sensiti vity analyses. 



it has been shown that the sampling design is one of the key issues (*S altelli et al 



I2OOO; Fang et al., 2006). Moreover, to avoid the problem of huge calculation time, 
it is often useful to replace the computer code by a ma t hematical approximation , 
called a su rrogate model or a metamodel ("Simp son et al. , 200 Ibt Fang et al. . l2006t 
Volkova et al.. 2008). The optimal exploration of the variation domain of the input 



variables is therefore especially important in order to avoid non-informative sim- 
ulation points ("Sobol', ' 1976'; 'Simpson et al.l l2001at iBursztvn and Steinber 2 12006I: 
[fooss et al., 2010; Levy and Steinberg, 2010J. 



Thirty years ago. iMcKav et al. (19791) have introduced the concept of Latin hy- 
percube sampling (LHS) for numerical experiments. Compared to simple random 
sampling (SRS) which insured independence between samples, LHS ensures the full 
coverage of the range of the input variables. More precisely, LHS allows to accu- 
rately reproduce the one-dimensional projections of the input sampling design. In 
terms of uncertainty and sensitivity analyses, it has been theoretically and experi- 
menta l ly proved tha t LHS is more precise and robust than simple rand om sample 
dSteinl Il987t lOwenL 119921 ISaltelli et al.L I2OOOI; iHelton and DavisL l2003h . Moreover, 
in the last twenty years, several improvements have been propo s ed in order to op- 



Pistone and Vicariol 2010t Jourdan and Franco , 2010l) . 



^F^i ^1—1 

timize the space filling properties of LHS designs (IParkL I1994t IFang et al.L I2006t 



Our starting point is that the initial LHS algorithm supposes independence be- 
tween input variables while in some situations, this assumption is irrelevant. First, 
let us recall the two forms of dependencies between variables: the statistical and the 
physical ones. 

- The correlation coefficient is the simplest measure of the statistical dependence 

between two variables. For the SRS, the so-called joint normal transfor m method 

consi sts in inducing a correlation structure on the transformed marginals dKurowicka and Cooke , 
2006h . The rank correlation coefficient, based on the rank transformation (which 



turns each variable value to it s ran k in the sample), is known as a more robust 
measure. Ilman and Conoverld 19821) have introduced an algorithm to consider rank 
correlations between variables in LHS. Some limitations of these two dependence 
measures have led to the introduction of other statistical dependence modeling 
(copulae, vines, etc., see lKurowicka and Cookd d2006l) ). 
- Physical dependencies between variables can arise when a variable has a formal 
relation in f unction of other variables. Such input constraints have been studied by 
Borgonovo (l2008i) which has proposed a novel way to solve the sensitivity anal- 



ysis problem in presence of equality constraints. Another currently encountered 
physical dependence, which is the subject of this paper, concerns the existence 
of inequality relations between the variables. It is the case when one variable is 
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physically constrained to be larger (respectively smaller) than another. For exam- 
ple, a geometric parameter (radius, height, etc.) of two physical objects can be 
subject to a rigorous increasing order if one object is included inside the other. 

When building the sampling design, the inequality constraints have to be honored 
in order to avoid some physical incoherence in the input sets that will be run with the 
computer model. A first solution could be to sequentially simulate the input variables, 
allowing to bound one variable by another in order to enforce an inequality constraint. 
However, as we will see, this procedure affects the one-dimensional projections of the 
sampling design. We need a procedure which separates the effect of dependence (the 
inequality constraints) from the effects of marginal distributions (i.e. the probability 
laws defined for each variable). To attain this objective, we propose an algorithm 
which builds a LHS satisfying the inequality constraints. 

This paper is devoted to the detailed presentation of this algorithm, called the 
constrained LHS (cLHS). In the next section, we introduce this algorithm by giving 
some examples. We compare it with a SRS-based algorithm and illustrate the algo- 
rithmic performances. In the third section, we explain in details the cLHS algorithm. 
As the inequality constraints can be too stringent to find a cLHS, we derive a nec- 
essary and sufficient condition proving its existence from an initial LHS. Then, our 
methodology is applied on a real problem involving welding simulation models. A 
conclusion gives finally some prospects to improve the cLHS algorithm. 



2 The sampling techniques 



The goal of the sampling step is to generate a matrix X" = {x^f)i=i..nj=i..p, where 
n is the number of experiments and p is the number of variables. The most common 
sampling met hod is indisput ably the pure Monte Carlo (i.e. SRS), mainly because of 
its simplicity dGentlel l2003h . It consists of randomly sampling n independent input 
variables. However, it is known to have poor space filling properties: SRS leaves 
large unsampled regions and can propose too close points. An example of a SRS is 
presented on Figure[T](a). 



2. 1 Latin hypercube sampling 



McKay et al. ( 19791) suggested an alternative method of generating X" that they called 



Latin hypercube sampling (LHS) which is an extension of stratified sampling. LHS 
ensures that each of the input variables has all of its range represented. Let the range 
of each variable Xj, j = 1 ... p, be simultaneously partitioned into n equally probable 
intervals. We note X" the n-sample of the variable Xj. A LHS of size n is obtained 
from a random selection of n values — one per stratum — for each Xj. Thus we 
obtain p n-tuples that form the p columns of the n x p matrix of experiments X" 
generated by LHS: the i''' line of this matrix contains the p input variables and will 
correspond to the code execution. Once a point is selected in an interval, no other 
point could be selected in this interval (see Figure [T](b)). Let us remark that the parti- 
tion into equally probable intervals allows to take into account non uniform densities 
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of probability like a normal distribution for example. Figure [T] shows 10 samples 
of two random variables obtained with SRS and LHS schemes. We can see that the 
result of LHS is more spread out and does not display the clustering effects found in 
SRS. 

Mathematically, if Xi, . . ., Xp are mutually independent random variables with 
invertible continuous distribution functions Fj, j ~ l,...,p, respectively, then the 
LHS i-th sample for the j-th variable can be created as 



J') 



,('■) 



'■) 



(1) 



where the Ttj are independent uniform random permutations of the integers {1,2,. ..,«}, 
and the ^j'^ are independent [0, 1] random numbers independent of the TZj. 



(a) Simple Random Sampling 



(b) Latin Hypercube Sampling 



0.0 0.2 0.4 0.6 0.1 
X1 



1.0 0.0 0.2 0.4 0.6 0.8 1.0 

XI 



Fig. 1 Examples of two ways to generate a sample of size fi = 10 from two variables X = [Xi ,X2] where 
Xi has a uniform distribution [0, 1] and X2 has a normal distribution ,_/V (0, 1). 



2.2 Constrained simple random sampling 



To take into account inequaUty constraints between variables, the simplest approach 
is based on SRS and consists on bounding one variable by another in order to enforce 
the inequality constraints. This approach is called the constrained Simple Random 
SampUng (cSRS). 

In Figure|2] we see the effect of an inequality constraint between two variables in 
terms of bivariate plots. Of course, the introduction of the truncation during the sim- 
ulation creates statistical dependences between X" and . This correlation depends 
on the distribution functions of Xi (called Fi) and X2 (called F2). In our example. 
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the correlation coefficient p{Xi,X2) is worth 31%. In Figure |2] the one-dimensional 
marginal projections of the samples are also shown. For the cSRS, the one-dimensional 
marginal of does not correspond to F2 anymore but to a transformed distribution 
F2 (which depends on Fi). 




(a) Simple random sampling (SRS) (b) Constrained simple random sampling (cSRS) 

Fig. 2 Comparison between Monte Carlo samples of size « = 100 from two variables X = [Xi ,X2] where 
Xi ~'2<[0, 1] andA:2 ~ '^^[O,!]. The inequality constraint forcSRS isXi <X2. 



This problem becomes more dramatic when the input dimension increases and 
when several sequential inequality constraints have to be satisfied, as for example if 
Xi < Xj+i for i — I,... ,p~l and p is large. Figure |3] (a) shows an example of cSRS 
of p = 10 variables with such increasing constraints. In this graph, each observation 
is represented as a line. Because of the sequential algorithm starting at Xi , all the 
curves are concentrated near the upper bound curve of the variables. Figure [3] (b) 
clearly reveals that the sampling of the first variable Xi is adequate with its uniform 
distribution, and that the samples of the following variables (X2 to Xiq) progressively 
take place in the upper region of their variation range. 

In summary, in some practical situations, users would like to simulate samples 
which follow all one-dimensional marginals and which take into account some in- 
equality constraints. The following section proposes such an algorithm. 



2.3 Constrained Latin hypercube sampling 



In order to follow all one-dim ensional margin als, our sampling procedure uses LHS. 
Our method, first proposed in IPeteleti (120071) . consists in doing permutations on an 
initial LHS to enforce the desired monotonic constraint. It is based on the fact that 
permuti ng two values of a variabl e in a LHS does not break the LHS structure of the 
sample dlman and Conoven,ll982l) . An appropriate algorithm scans the starting LHS 
to find the couples of values that violate the monotonic constraint. Then the algorithm 
finds and executes the combinations of permutations which have to be done to satisfy 
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the inequality constraint between X, and for the n experiments. Details of the 
algorithm are given in section [3] 

Figure|4]shows the work done on a couple of parameters on which an increasing 
constraint is enforced. The distribution of X\ and X2 are kept uniform (the slight 
variation of the height of one class in the histogram of X2 is not significant). In the 
bivariate plots, we see that the cLHS constraints (the increasing inequality and the 
honoring of all one-dimensional marginals) tend to gather the sample points along the 
inequality frontier line X\ ~X2. It appears that the severity of the inequality constraint 
effects strongly depend on the one-dimensional marginal distributions Fi and F2 of 
the variables. The limit case is illustrated on Figure |4](b) where Xi and X2 have the 
same one-dimensional marginal distributions Fi and F2 (then the same upper and 
lower bounds). In such a case where half the area of the bivariate plot is forbidden all 
the points are located on this frontier Une. This effect result s of a too seve re constraint 
and reveals the need of a constraint intensity measurement. I Peteletl(l2007l) has defined 



this constraint intensity measurement as the ratio between the triangular forbidden 
area (St) in the bivariate plot and the rectangular area (Sr) of the domain defined by 
the upper and lower bounds of the variables: 

The 7 measurement can be used if each of the input variable has some upper 
and lower bounds, i.e. if the support of their distribution function is defined on a 
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(a) 



(b) 



Fig. 4 Comparisons between constrained Latin liypercube samples of size n = 100 from two variables 
X=[Xi,X2] withXi ^^/[0,l] and the inequality constraint Xi < X2: (a) X2 ~ [0,2]; (b) X2 ~ "^^lO, 1.1]. 



bounded domain. In the general case of an inequality constraint between X, and Xj, 
respectively defined on [bi,hi] and [bj,hj], we obtain 



ST{Xi,Xj) 
SR{Xi,Xj) 



— for the constraint X- < Xi , 

2{h,~b,){hj~bj) ^' 

(h ■ - bif ^ ' 

— , ' 7 for the constraint Xi > Xi . 

2{hi^b,){hj~bj) ^ 



The intensity constraint measurements for the Figure |4]cases are worth y(Xi , X2) = 
25% for (a) and 7(^1, X2) = 45.5% for (b). With some heuristic arguments. IPetelel 
(I2007I) has found a nearly linear link between 'Y{Xj,Xj) and the correlation coefficient 

p(X,^/)for 7(X,Xy) e [0,0.3]: 



piXi,Xj)^2.nS 7{Xi,Xj) 



(4) 



As 7 is positive, the correlation coefficient will be always positive. For example, 
this relation shows that if the inequality constraint is kept smaller than 15%, the 
correlation between the variables will be smaller than 40%. 

For the same cases than in Figure |3] Figure |5] shows the cLHS of p = 10 vari- 
ables with sequential increasing constraints. As before, in this graph, each observa- 
tion is represented as a line. At present, the curves correctly fill the variation ranges 
of all the variables. The constraint intensity measurements for Figure|5](a) are worth 
Y{Xi,Xj^i) — 18.75% for / = 1, ... ,9. This value is rather suitable: correlations be- 
tween variables are smaller than 52% (value obtained thanks to Eq. (|4|). For Fig- 
ure|5](a), the constraint intensity measurements increase from y{Xi ,^2) = 33.33% to 
7(X9,Xio) =45.45%. 

When the upper and lower bounds of the variables are similar, the cLHS tends to 
give homothetic translated trajectories, as shown by Figure |5](b). In the limit case. 
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if all the X,'s have the same one-dimensional marginal distributions (then the same 
upper and lower bounds), the obtained curves are parallel and regularly spaced be- 
tween the lower bound curve and the upper bound curve. This is one of the drawback 
of our algorithm, caused by the imposed LHS structure. Moreover, the feasability of 
the cLHS depends on the bound values of the constrained variables. For example, for 
the constraint Xj < Xj, the algorithm does not work if hi > hj or if bi > bj. From Eq. 
(|3]l, this implies that the constraint intensity measurement 7(X, ,X, ) is upper bounded 
by 0.5. 




(a) (b) 

Fig. 5 Constrained Latin hypercube samples of size n= 10 from p = 10 variables Xj, i = 1, . . . , 10, with 
Xi < Xi-^- 1 for (' = 1 , . . . , 9. The upper and the lower curves represent the bounds of the variation range for 
these 10 variables, (a) X, ~ ^/[0+ ^,2+^]; (b) X, ~ f/[0, 1 + i^]. 



3 The constrained Latin hypercube sampling algorithm 

In the following, we explain the cLHS algorithm in the case of a strict increasing 
constraint between two variables Xi and X2 (with distribution functions Fi and F2 re- 
spectively). The developments for the strict decreasing constraint case are exactly the 
same, by inverting the inequalities sense. Moreover, the extension of our algorithm to 
non strict inequality constraints is straightforward. 

For the increasing constraint case, we assume the following hypotheses: 

- X = {Xi,X2) is defined on a bounded domain ^ G M?. The support of Fj for 

{1,2} is [bj.hj]; 

- The bounds are subject to the following inequalities: 

bi < b2 and hi < /12 . (5) 

These inequalities seem natural: if we impose some increasing constraint between 
Xi and X2, we hope that the same increasing constraints exist for their minimal 
and maximal bounds. 



9 



Let us define the matrix C" = C(Xj',X^) of size n x n: 

^" = i^ij)i=l..n,j=l..n~ : '■• : 



(6) 



where lj;>y = 1 if x > y and Ijoy = otherwise. C" is called the compatibility matrix 
between X" and X". This matrix allows to identify which combinations of elements 
of X" and X2 are incompatible, i.e. those with a decreasing relation. Therefore, the in- 
equaUty constraint between the two samples is honored if the diagonal of C" contains 
only I's. 

For our cLHS algorithm, if we choose to leave X" unchanged and give the pos- 
sibiUty to permute some elements of , we define our final objective as getting a 
sample Xj" such that 

L4=«, (7) 

!=1 

the compatibility matrix between X" and Xj". 

' i=\..n,J=l..n 

At present, it would be convenient to know if this objective can be achieved for a 
specific sample X". Let us define the sample vector 

S" = (50i=i..„= ( £c,v) . (8) 

Si gives the number of elements of XJ* which satisfy the constraint with ' . We 
also define {Si)^^^ „ the ordered sample of (>S',),=i..„. To insure that we can obtain a 
sample Xj" (by permutations of the X" elements) satisfying the increasing constraint, 
the following assertions have to be true: 

- The smallest element of Xj have one or more than one smaller elements in X", 
which is equivalent to say that nain(X2 ) > min(X"), then to say that >Si > 1; 

- The i*^- smallest element of Xj have i or more than i smaller elements in X", 
which is equivalent to say that Si > i; 

- The n^^-smallest element of X^ have n or more than n smaller elements in X", 
which is equivalent to say that 5„ > n; 

From these assertions, we obtain the following result: 

Proposition: If (X^Xj ) is a LHS, the inequality 

min[(5i ••• •••«)] >0 (9) 

is a necessary and sufficient condition to guarantee the existence of a sample Xj" 
such that (XfjXj") satisfies the increasing constraint Xi < X2, where Xj" is obtained 
from permutations of the elements o/Xj . 
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Therefore, the first step of our methodology will be to test this criterion. If equa- 
tion (|9|l is not verified for a chosen X" — (Xj^Xj ), a new sample for X2 (keeping 
the LHS property for X") is created and the scan starts again. Our hypotheses on 
the bounds of Xi and X2 (Eq. (|5])) guarantee that a LHS satisfying the increasing 
constraint exists. 

We start from an initial LHS X" = (X",X2 ), with a compatibility matrix C, satis- 
fying the existence criterion We want to obtain the LHS X" = (ZfjXj") (with a 
compatibility matrix C') satisfying the increasing constraint between Xi and Xj. Our 
objective is therefore to obtain the result of equation Q. 

Let us note X" the reverse ordered sample of X". This vector X" contains the 

elements Xj^' > xf'^ >■■ ■> i j"' . We put in the sample vector A" a sequence of indices: 
the indices in X" of the X" elements. Mathematically, it follows that 



Our permutation algorithm is based on the treatment of the X" elements in a 
sequential manner (because the constraint is more difficult to be satisfied by the first 
values of X"). We describe the algorithm with the following four steps. 

Algorithm cLHS: 

1. Initialisation: C'" = C" and B' = (1 • • • n). 

2. For i — I ...n: 

- We put in the vector B the indices in of the elements compatible with x\ : 



- We randomly choose an element in B and put it in B'^. 

- The index B'^ corresponds to the one that will be permuted in . We turn to 



3. The vector B' contains the indices that will be used to make the permutations in 
^2 . We obtain the new sample of the vaiiable X2: 



k=l 



For j = I ...n: 

if c'j^. = 1 then B^ = k and k = k+ 1 




(in order to block up the index 



For / = 1 . . . n : c',. = 




(11) 



4. Finally, the permutation matrix C " is calculated with X" and Xj" by equation ^ 
in order to test the equality of equation (|7]i. 
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End of algorithm 

The extension of the CLHS algorithm to the multivariate case X = {X\ ,Xp) 
with Xj < Xj^i for j = l,...,p ~ 1 is straightforward and is done in a sequential 
manner. We first simulate a LHS X" = {X", ...,Xp). Leaving X" unchanged, we se- 
quentially build with the cLHS algorithm x'" from x'"_^ and X" for j ^2,...,p. At 
each step /', before applying the algorithm, the criterion ^ is tested. If this criterion 
is not verified, a new LHS for X" is created, and so on until the criterion is verified. 



4 The cLHS algorithm in an example 



We propose a simple example with a sample of size n = 6 and two variables Xi and 
X2 subject to a decreasing constraint Xi > X2- Points are uniformly sampled on the 
domain ^ = [20,30] x [16,26] of (Xi,X2). We simulate an initial LHS (Figure|6](a)) 
and obtain the following design matrices: 



X'{ = 



/23.98\ 
26.91 
26.52 
2L99 
29.23 

V2I.IO/ 



X^^ 



then C" = 



/I 1 1 1 0\ 
111111 
1110 10 
111111 
111111 

VO 1 1 1 0/ 



/22.18\ 
20.45 
23.77 
18.31 
16.45 
^25.49/ 

As the diagonal of this matrix has a null term, the equation (|7]) is not fulfilled. Then, 
the cLHS algorithm has to be applied in order to obtain a LHS satisfying the con- 
straint. 

First, we test the existence criterion. We obtain S" = (46 46 6 3) and S" = 
(344666). The existence criterion (Eq. ^) is then fulfilled. 
Second, we apply the algorithm and obtain the following result: 

/20.45\ 
25.49 
22.18 
18.31 
23.77 
VI6.45/ 



Xo 



(12) 



X" has not been modified while elements of have been permuted to obtain Xj", 
which is a sample satisfying the decreasing constraint. Other x'2 samples could be 
found, the choice made during the cLHS algorithm being random. Figure|6](b) shows 
our final sampling result. 



5 An appUcation case: welding thermomechanical models 

The robust increase in computer power has tremendously contributed to a growing 
fad for welding simulation. The industrial requirements are more and more numer- 
ous: supports to develop new processes, control of mechanical welding effects (in 
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(a) Initial LHS (b) Final LHS 



Fig. 6 Illustration of the cLHS algorithm (n = 6) with a decreasing constraint between two uniformly- 
distributed random variables Xi and X2 (7(^1,^2) = 18%). The diagonal line corresponds to the frontier 
lineXi =X2. 



particular, residual stresses and distortions), argument in a nuclear safety analysis re- 
ports, etc. Thus, through the use of high-performance computers and advanced mod- 
els, numerical simulation is expected to become an important tool for innovation in 
welding engineering. 

However, running a welding simulation model requires a large number of inputs 
- about 500 - including for example meshing inputs, boundary and initial conditions 
as well as material properties and process parameters, and generates several outputs, 
including spatial distributions of displacements and residual stresses in the weldment. 
In particular, among inputs, the determination of material properties is one of the key 
problems of welding simulation. The features of material properties are that they are 
dependent on temperature and that their full characterization is very expensive, often 
difficult or even sometimes impossible. In this context, the global sensitivity analyses 
of the numerical welding simulation model allows to determine which material prop- 
erties are the most sensitive in a numeric al wel ding simula tion and in which range of 
temperature dPetelell 120071 fPetelet et all I2OO6I: IXsserin et all l2009l) . 

Let us show the application of our methodology on the range of steel material. 
Five input variables are the mechanical properties used by the model: Young's mod- 
ulus, thermal computation coefficient, Poisson's ratio, yield strength and hardening 
modulus. Because of their dependence on temperature, it has been required to sam- 
ple each material property at a discrete set of temperatures: 7 levels are chosen from 
20"C to 1100"C. We obtain 35 = 5 x 7 input variables, each following an uniform 
distribution defined by its minimal and maximal bounds (taken from the literature). 
Moreover, some material properties used in this model are monotonically decreas- 
ing as function of temperature. Therefore, the constrained Latin hypercube sampling 
strategy, described in this paper, can be used to generate the input design. This strat- 
egy allows our sampled variables to honor their uniform repartition defined by their 
minimal and maximal bounds, that is to say sampUng in the physical bounds. 
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To illustrate this application, we present the sampling of the Young's modulus. 
Figure |7] shows the result on this 7-dimensional variable which follows a decreas- 
ing constraint. These curves present three randomly selected materials among the 
800 created and the bounds of the domain. For the sen sitivity analysis proc ess, the 



Young's modulus is represented by only 7 parameters dAsserin et al. ■ l2009l) . How 



ever, one should keep in mind that for the mechanical computation, the curve rep- 
resents truly the considered dependence of this modulus because the algorithm uses 
intermediate values according to a piecewise linear interpolation. 



3.0E+11 
2.5E+11 
2.0E+11 
1.5E + 11 
1.0E+11 
5.0E+10 
O.OE+00 

100 200 300 400 500 600 700 800 900 1000 1100 

Fig. 7 Example of n = 3 samples of the Young's modulus obtained by cLHS (dashed lines). The upper 
and the lower curves (solid lines) represent the bounds of the domain for these p = l variables. 

It would be interesting to present the sensitivity analysis results obtained by em- 
ploying the classical LHS (unconstrained) and by employing the constrained LHS. 
However some model calculations do not converge when using non physical evolu- 
tion of the material properties (i.e. unconstrained LHS). Moreover, performing model 
calculations with non monotonic evolutions would be an aberration from a physical 
point of view. 




6 Conclusion 

In this paper, we have proposed a new algorithm (called cLHS) allowing to obtain 
samples of several variables constrained by some inequality relations. This situation 
can frequently arise in application cases while very few works have been devoted to 
this issue. The cLHS algorithm allows to satisfy the inequality constraint while leav- 
ing unchanged the one-dimensional marginals that we have defined for each variable. 
In order to honor these one-dimensional marginals, the LHS-based technique has 
been preferred. To our knowledge, this inequality constraint problem has not been 
studied for the LHS building issue. 

We have shown the interest of this algorithm in an application case involving 
welding simulation model. The cLHS algorithm has proven its efficiency to sample a 
multi-dimensional variable taking under consideration its physical nature. 

The current cLHS algorithm has one main drawback. When the minimal (respec- 
tively maximal) bounds of the two constrained variables move closer, the space filling 
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properties of the sample points deteriorate: the sample points are gathered along the 
inequality frontier line. The cLHS is therefore efficient if the bounds between the 
variables are sufficiently distant. A constraint intensity measurement, noticed 7, has 
been defined in order to quantify this effect. Moreover, a linear relation has been 
proposed between 7 and the correlation coefficient of the constrained variables. This 
allows to a priori know (from the variable bound values), the effects of the inequality 
constraint in terms of correlation of the simulated sample. 

In a future work, it will be interesting to quantify this phenomenon by linking 
a space filling measure (as the discrepancy) with the bound values. More generally, 
if we want to suppress this undesired effect, the LHS framework has to be left out. 
A first idea would be to work with entropy-based designs by optimizing the entropy 
on one-dimensional marginals. Developing design optimization algorithms under in- 
equality constraints would be an interesting research way. 

If more than two variables are under study, the cLHS algorithm is limited to in- 
equality constraints defined in a sequential order A complex inequality constraint 
could involve more than two variables (e.g. Xi < X2 +X3) or could be non sequential 
(e.g. Xi < X2 and Xi < X3). If such situations are identified in specific application 
cases, there is no doubt that some extensions of the cLHS algorithm are possible. 
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